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All evolutional paths for all admissible initial conditions of FRW cosmological models with dissi- 
pative dust fluid (described by dark matter, baryonic matter and dark energy) are analyzed using 
dynamical system approach. With that approach, one is able to see how generic the class of so- 
lutions leading to the desired property - acceleration - is. The theory of dynamical systems also 
offers a possibility of investigating all possible solutions and their stability with tools of Newtonian 
mechanics of a particle moving in a 1-dimensional potential which is parameterized by the cosmo- 
logical scale factor. We demonstrate that flat cosmology with bulk viscosity can be treated as a 
conservative system with a potential function of the Chaplygin gas type. We characterize the class of 
dark energy models that admit late time de Sitter attractor solution in terms of the potential func- 
tion of corresponding conservative system. We argue that inclusion of dissipation effects makes the 
model more realistic because of its structural stability. We also confront viscous models with SNIa 
observations. The best fitted models are obtained by minimizing the function which is illustrated 
by residuals and levels in the space of model independent parameters. The general conclusion 
is that SNIa data supports the viscous model without the cosmological constant. The obtained 
values of statistic are comparable for both the viscous model and ACDM model. The Bayesian 
information criteria are used to compare the models with different power law parameterization of 
viscous effects. Our result of this analysis shows that SNIa data supports viscous cosmology more 
than the ACDM model if the coefficient in viscosity parameterization is fixed. The Bayes factor is 
also used to obtain the posterior probability of the model. 

PACS numbers: 98.80.Bp, 98.80.Cq 



I. INTRODUCTION 

Several astronomical observations, since observations of distant supernova type la P, Q and cosmic microwave 
background anisotropy measurements [3, 0] , have indicated that our Universe is currently undergoing an accelerating 
phase of expansion. If we assume that the Friedmann equation with perfect fluid source is valid then our Universe 
is in acceleration phase due to the presence of matter with negative pressure violating the strong energy condition 
(dark energy) . While there are several candidates for description of dark energy like an evolving scalar field (referred 
to as quintessence) Q , generalized Chaplygin gas [1, 0] or phantom dark energy |1, [p | the most obvious candidate 
for such a component of the dark energy seems to be the cosmological constant [I^]- Hence we obtain the most 
simple explanation of current universe which require two dark components, namely non-relativistic dust (dark and 
baryonic matter) contributing about 30% of the total energy density (which is clustering gravitationally at small 
scales) and second component (dark energy) which dominate at large scales and is smoothly distributed with constant 
density. Unfortunately it is only some kind of effective theory which offers description of present observations rater 
then its understanding. For example it remains to understand why the value of cosmological constant obtained 
from observations of SNIa is so small in comparison with the vacuum energy (Planck mass scale). It is the main 
reason for looking for alternative explanation of observational data of distant supernovae. There have appeared many 
models that make use of new physics which may arise as a consequence of embedding of our universe in a more 
dimensional universe which accelerates due to some additional term in the Friedmann equation rather than dark 
energy contribution pT| . 

In this work we discuss the possibility of the dark energy being represented by viscous fluid which is characterized 
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by the presence of bulk viscosity coefficient ^ = —{l/3)dpcs/dH in the effective equation of state. Because we assume 
homogeneity and isotropy, only bulk viscosity effects may be important under such symmetry condition. FRW models 
filled with a non-ideal (dissipative) fluid were investigated in several papers in which bulk viscosity was parameterized 
by an additional term in pressure p — > pcff = p — 3i^i?, where ^ is the bulk viscosity coefficient, usually parameterized, 
following Belinskii, as power law ^(p) = ap™, where p is the energy density and a and m are constants jl^ ITsl [l^ ITsI . 
While we know since the classical papers of Murphy and Heller that bulk viscosity effects can generate an 
accelerating phase of expansion the idea that the present acceleration epoch is driven by some kind viscous fluid is 
being explored recently 18]. Wilson et al. [l^ discussed a cosmology in this cold dark matter particle decay into 
relativistic particles and argued such a decay could lead naturally to bulk viscosity modelled by effective negative 
pressure. 

The main aim of the paper is to investigate the possibility of explaining the present accelerating expansion of 
the Universe in terms of bulk viscosity with the background of models with Robertson- Walker symmetry, basing on 
analogy with the cosmology with the Chaplygin gas and dynamical system methods. We show that the non-flat FRW 
dissipative models can treated as a perturbed dynamical systems of the Newtonian type. 



The Chaplygin gas conception is taken from aerodynamics. Kamenshchik at al. [20] noted that from this equation 
of state can interpolate two stages - the stage of matter domination and the stage of dark energy domination. Then it 
was generalized to the form of the generalized Chaplygin gas p — —A/p" by taking a different from one. The reason 
was to increase the degree of generalization without any physical justification. 

To have the physical meaning of this class of models we proposed to consider the FRW model with bulk viscosity 
which is formally indistinguishable from the flat FRW model with the generalized Chaplygin gas. The advantage of 
the FRW model with bulk viscosity is that its parameter possesses the physical interpretation - dissipation in the 
sense of bulk viscosity although the integral of energy is conserved and this model belongs to the class of models 
called conservative models described in terms of the potential function. 

One can distinguish between conservative and dissipative dynamical systems (divergence of the vector fleld is zero 
for conservative, larger then zero for dissipative). We show that the flat FRW with viscous fluid is conservative in the 
above sense. If we consider non flat models then the analogy between the flat FRW model with bulk viscosity and 
the Chaplygin gas model vanishes. 

In investigation of dynamical effect of bulk viscosity we use methods of dynamical system theory which offers the 
possibility of investigating all evolutional paths for all admissible initial conditions. In this approach the key point is 
the visualization (geometrization) of dynamics with the help of phase space. We consider dissipative cosmology (as 
well as conservative) in terms of a single potential function which determines the dynamics on a 2-dimensional phase 
space (a, a), where a is the scale factor, '= d/dt, and t is the cosmological time. 

Inspired by the fact that viscous fluid possesses negative pressure and in flat models can be modelled as a Chaplygin 
gas we have undertaken the simple task of studying the FRW viscous cosmology. The dissipative cosmological model 
have at least a few significant advantages from the theoretical as well as observational point of view: 

1. they can describe a smooth transition from a decelerated - matter dominated expansion of the Universe to the 
present epoch of cosmic acceleration; 

2. they offer the possibility of unification on phenomenological ground both dark energy and dark matter; 

3. they can be treated as a natural extension of the CDM models in which effects of dissipation are unified; moreover 
from confrontation with SNIa data they are a serious alternative to the concordance ACDM model; 

4. in the flat cosmology they can be identified with conservative models which are structurally stable, i.e. such 
that a small perturbation does not change their global scenario. 

Organization of the text is the following. In section |TT] we describe all classes of FRW conservative models with 
dark energy in terms of potential function for fictitious particle moving on the constant energy level. In section IIIII 
the dynamics of conservative dark energy models with dissipation is investigated in tools of qualitative theory of 
differential equation in the 2-dimensional phase space. The full analysis of dynamics of viscous models in terms of 
perturbed conservative systems is performed for the constant bulk viscosity coefficient. We also define the distance the 
conservative and perturbative systems using the Sobolev metric. Section |TV] contains some observational constraints 
from distant supernovae type la. These observational constraints indicate that unified model for dark energy and 
dark matter through the employment of dissipation is favored over the ACDM model by the Bayesian information 
criteria of model selection. 
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TABLE I: Different dark energy models in terms of the potential function parameterized by the scale factor or redshift z 
(l + z^a-i). 



model 


potential function 


independent parameters 
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PhCDM model 
(phantom wx ~ ~4/3) 


V{a{z)) = -i{!^m,o(l + z)+ Qphfiil + zy'-^ + fifc.o} 




oivv^uivi moQei 
bouncing with A 
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^ ^ o o \ n- 


Randall- Sundrum 
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linearized at 2: = 
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models with generalized 
Chaplygin gas 

px = --^ 

Px 

and baryonic matter 


1 

V{a{z)) = -i{nb,o + nch,o(l + 2)-" + (1 - + 2)=*''+") 
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II. COSMOLOGICAL MODELS WITH DARK ENERGY AS A CONSERVATIVE SYSTEM 

The FRW dynamics of a broad class of dark energy models can be represented by the motion of a fictitious particle 
in one dimensional potential V{a), where a is the scale factor. The equation of motion assumes the form of a 
2-dimensional dynamical system of a Newtonian type 

dV 
oa 

where V{a) = and effective energy density satisfies the conservation condition 

Pcff = -3i?(pcff + Pcff), (2) 

where H ~ a/ a is the Hubble function and a is the scale factor expressed in units of its present value ao = 1. Therefore 
different dark energy universe models filled with dust matter and dark energy satisfy the general form of the equation 
of state 



Px = wx{a)px, (3) 



where 

PcS = 0+Px, 



Pcff = Pm^oa ^ + px,oa ^ exp 



are the effective pressure and energy density, respectively. The index zero denotes that the quantities are calculated 
at the present epoch. They can be classified in terms of the potential function (see Table |T|. In this table we used 
redshift z instead of the scale factor a because of the relation 1 + z — and dimensionless density parameters for 
each component U^ q — pijiH^. 

There are different advantages of representing the evolution of the Universe in terms of a 2-dimensional dynamical 
system of a Newtonian type. Firstly, we obtain unified description of the broader class of cosmological models with 
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a general form of the equation of state. Secondly, the shape of potential which determines the critical points and 
their stability can be reconstructed immediately from the SN la data due to simple relation between the luminosity 
distance (ii(z) and the Hubble function in the case of the flat model. Thirdly, the potential function can be useful 
as an instrument of probing of dark energy models from observations because as opposed to the coefficient of state 
wxiz) it is related to the luminosity distance dL{z) by a simple integral. 

For our further analysis it is important that all these systems are conservative ones. They are given on zero energy 
level which is conserved 

'^ = IpI + ^ 0, (4) 

where pa = a. 

The dynamic of the models under considerations is reduced to the dynamical system of the Newtonian form 



dV^ . (5) 

where x = a/ao and the above system has the first integral (energy) in the form 

yl + V{x)^0. (6) 

It is important that the critical points of ([5]) as well as their character can be determined from the geometry of the 
diagram of the potential function only. 

The system ([SJ has critical points localized on x-axis because x = xq and {dV/dx)\xo = at the critical point. 
Therefore they are representing the static solution. The domain admissible for classical motion is Dq = {x- V{x) < O}. 

For the conservative system ([S]) it is useful to develop methods of qualitative investigations of differential equation 
[2l| . The main aim of this approach is to construct phase portraits of the system which contain global information 
about the dynamics. The phase space {x,y) offers a possibility of natural geometrization of the dynamical behavior. 
It is in simple 2-dimcnsional case structurized by critical points or non-point closed trajectories (limit cycles) and 
trajectories joining them. Two phase portraits are equivalent modulo homeomorphism preserving orientation of the 
phase curves (or phase trajectories). From the physical point of view critical points (and limit cycles) are representing 
asymptotic states of the system or equilibria. Equivalently, one can look at the phase flow as a vector field 



dV 
ox 



(7) 



whose integral curves are the phase curves. 

Thanks to Andronov and Pontryagin [22] , the important idea of structural stability was introduced into the multi- 
verse of all dynamical systems. A vector field, say /, is a structurally stable vector field if there is an e > such that 
for all other vector fields g, which are close to / (in some metric sense) ||/ — g\\ < s, f and g are topologically equiv- 
alent. The notion of structural stability is mathematical formalization of intuition that physically realistic models in 
applications should posse some kind of stability, therefore small changes of the r.h.s. of the system (i.e., vector field) 
does not disturb the phase portrait. For example motion of pendulum is structurally unstable because small changes 
of vector field (constructed from r.h.s. of the system) of friction type g = [y, —dV/dx + ky]^ dramatically changes the 
structure of phase curves. While for the pendulum, the phase curves are closed trajectories around the center, in the 
case of a pendulum with friction (constant) they are open spirals converging at the equilibrium after infinite time. We 
claim that the pendulum without friction is structurally unstable. Many dynamicists believe that realistic models of 
physical processes should be typical (generic) because we always try to convey the features of typical garden variety 
of the dynamical system. The exceptional (non-generic) cases are treated in principle as less important because they 
interrupt discussion and do not arise very often in applications 23]. 

In the 2-dimensional case, the famous Peixoto theorem gives the characterization of the structurally stable vector 
field on a compact two dimensional manifold. They are generic and form open and dense subsets in the multiverse of 
all dynamical systems on the plane. If a vector field / is not structurally stable it belongs to the bifurcation set. 

The space of all conservative dark energy models can be equipped with the structure of Banach space with the 
metric. 

Let Vi and V2 be two dark energy models. Then distance between them in the multiverse is 

d(yi,V2) = max] sup {Vi^oi: - V2,x\, sup \Vi^rx - V2,xx\\, (8) 
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FIG. 1: Three types of behavior in the neighborhood of a critical point: a saddle, a center and a cusp. 

where £^ is a closed subset of configuration space. Of course, multiverse of all dynamical systems of a Newtonian type 
on the plane is infinite dimensional functional space and the introduced metric is a so-called Sobolev metric. 

While there is no counterpart of the Peixoto theorem in higher dimensions it can be easy to test whether planar 
polynomial systems, like in considered case, have structurally stable phase portraits. For this aim the analysis of 
behavior of trajectories at infinity should be performed. One can simply do that using the tools of Poincare 
construction, namely by projection trajectories from the center of the unit sphere = {{X, Y, Z) G R'^ : + Y^ + 
= 1} onto the {x,y) plane tangent to 5^ at either the north or south pole. 

The vector field / is structurally unstable if 

1. there are non- hyperbolic critical points on the phase portrait, 

2. there is a trajectory connecting saddles on the equator of S^. 

In opposite cases if additionally the number of critical point and limit cycles is finite, / is structurally stable on S'^. 

Let us consider a 2-dimensional dynamical system of a Newtonian type ([5]) . There are three cases of behavior of 
the system admissible in the neighborhood of the critical point {xq, 0) : — dV/dx\xo = 

• If {xq, 0) is a strict local maximum of V{x), it is a saddle point; 

• If {xq, 0) is a strict local minimum of V{x), it is a center; 

• If {xq, 0) is a horizontal infiection point of V{x), it is a cusp. 

All these cases are illustrated in Fig. [1] 

It is a simple consequence of the fact that characteristic equation for the linearization matrix at the critical point 
{tT A = 0) is + detA — 0, where det^ = d'^V/dx'^\xa- Therefore the eigenvalues are real of opposite sign for saddle 
point, and for centers - if they are non-hyperbolic critical points - purely imaginary and conjugated. 

Finally, if the set of all conservative dark energy models has potential function V{x) such that the number of static 
critical points and cycles is finite and there are no trajectories connecting saddles then / € C^{A4), r > 1 contains 
open, dense subset of all smooth vector fields of class C" on the phase space. 

In Fig.[2]the phase portrait for the ACDM model is presented on a compactified projective plane by circle at infinity. 
Of course it is structurally stable. Therefore, following Peixoto theorem, it is generic in the multiverse M of all dark 
energy models with 2D phase space because such systems form open and dense subsets. 
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FIG. 2: The phase portrait for the ACDM modeL On the phase portrait we have a single saddle point in the finite domain and 
four critical points located on the circle at infinity. They represent initial singularity (a; = 0, i = oo) or de Sitter Universe deS. 
The trajectory of the flat model (fc — 0) divides all models in to two disjoint classes: closed and open. The trajectories situated 
in the region II confined by the upper branch of the fc = trajectory and by the separatrix going to the stable de Sitter 
node deS+ and by the separatrix going from the initial singularity to the saddle point S correspond to the closed expanding 
universe. The trajectories located in the regions I and II (which corresponds to the open universes) are called infiectional. 
Quite similarly, the trajectories situated in the region III correspond to the closed universes contracting from the unstable de 
Sitter node towards the stable de Sitter node. The trajectories running in this region describe the closed bouncing universes. 



It is also interesting that the phase space of dark energy model can be reconstructed from SN la data set (Gold 
Riess sample) and is topologically equivalent to the ACDM model. Fig. [3] represents the potential function 



d dUz) 
dz 1 + z 



(9) 



reconstructed from dL{z) relation - luminosity distance rf^ as a function of redshift z : 1 + z = a ^. Such a 
reconstruction is possible due to the existence of the universal formula 

dLiz) _ f dz' ^^^^ 



1 + z J H{z'y 
for the flat model. 

In principle it is good news for the ACDM model and others which are close to it in the multiverse in the metric 
sense. 

On the phase portrait the phase space of the Cardassian model is presented (Fig. |4|) . In the case of single fluid, the 
presence of additional terms in the modified Friedmann equation is equivalent to extra contribution to the effective 
energy density in the standard Friedmann equation. On the phase portrait there is a degenerated point x = oo, y = 
located at the infinity. Let us note that there is no de Sitter attractor at late times. This obstacle was addressed in 
Ref. [H,!!^. 

Analogous situation appears if we consider negative value of the exponent n in the Cardassian models. It is well 
known that if we consider dust matter only Cardassian models with n < are equivalent to phantom models. The 
phase portrait for the phantom model {px — — (4/3)px) with dust is illustrated in Fig. [51 In both cases, Fig. 2] and 
Fig. [5l the systems are structurally unstable due to the presence of degenerate critical points at infinity. The models 
with de Sitter phase of evolution at late time are distinguished in the multiverse of dark energy models. 

Our main objection addressed to the Cardassian models is that they do not offer description of dark energy epoch 
in terms of a structurally stable model. Cardassian cosmologies belong to the non-generic, in our terminology, class 
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FIG. 3: The potential function as a function of the scale factor expressed in its present value ao = 1 for the reconstructed best 
fit model is given by the solid line. The confidence regions la and 2a are drawn around it. The phase portrait obtained from 
this potential (best fit) is equivalent to the ACDM model (see Fig. ^ . The value of the redshift transition estimated from SNIa 
data (the Gold Riess sample) is about 0.38 (see [53]). 




FIG. 4: The phase portraits for the Cardassian model n = —1 (left) and n = 2 (right) (see Table |l]). The Cardassian model 
with n = —1/3 corresponds to the phantom model with wx = const < — 1 (generally if we consider dust model equation of 
state Pcard = {n — l)pcard rcproduces the Cardassian term). 
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FIG. 5: The phase portrait for the phantom model. While at finite domain we find only a saddle point like for the ACDM 
model, at infinity, the trajectory is reaching the degenerate point a; = 0, i = oo. Note that this degeneration can be removed 
be redefinition of the phase variables and reparameterization of time. Then we obtain other system which is topologically 
equivalent to the ACDM model (see Fig. [2}. 



of fragile model of the present Universe. Similar models form the bifurcation set in the multiverse of dark energy 
models. 

There is also another example of the fragile model used in the cosmological applications ~ the bouncing model. 
Such a model recently appeared in the context of semi-classical description of the quantum evolution of the Universe 
in the framework of loop quantum gravity. In this scenario singularity separates classically allowed regions [27i] . The 
phase portrait of bouncing cosmological model is shown in Fig. [B] In such a cosmology one encounters a bouncing 
phase of evolution instead of big bang during which the potential function has a minimum. It is structurally unstable 
due to the presence of non-hyperbolic critical point on the phase portrait. In the next section we find some perturbed 
system which modifies this and makes it structurally stable. 



III. CLASS OF DISSIPATIVE DARK ENERGY MODELS 



Applying the analogy with pendulum with friction, we would like to incorporate dissipative processes into cosmo- 
logical model with dark energy. As a matter of fact, dissipative processes should be present in any realistic theory of 
the Universe 

The simplest way to include bulk viscosity effects is through Eckart's theory [29j, which postulates that the bulk 
viscous pressure is proportional to the expansion. Many authors used that theory to investigate bulk viscosity on the 
evolution of the Universe [H, [lB| . 

Basing on Israel and Stewart theory [sol . [sij of truncated viscosity, Belinskii [l^ studied effects of viscosity on 
cosmological models in terms of effective pressure 

PoS = p-3ap)H, (11) 

where p is the pressure of a fluid flUing the Universe, and 3^{p) = ap™ is the viscosity coefficient parameterized by 
energy density p = poff- 

In the context of viscous cosmology many authors used dynamical system methods in investigations of the stationary 
states and their stability [s^]. Several authors have studied cosmological models with viscosity [1^, which are also 
interesting in the context of inflation [33|. It is also interesting that viscosity may describe phenomenologically 
quantum effects [ssj . 
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FIG. 6; The phase portrait for bouncing (in the contemporary sense see Table cosmological model with A. On the s-axis, 
there are two critical points: center and saddle point. The global phase portrait (in contrast to the ACDM model) is structurally 
unstable because of the presence of a non-hyperbolic critical point - center. Note that all models (closed, open and flat) undergo 
a bounce, i.e., they begin their contraction from the unstable de Sitter node, reach a point of minimal contraction and begin 
expansion toward the stable node deS+. 



In this section the effects of dissipation are treated as a small perturbation of conservative dark energy models. We 
then investigate the influence of dissipation on global dynamics of fragile dark energy models. 
Let us start with the definition. 

Definition 1 A multiverse of the FRW dissipative models with dark energy is defined as a space of all 2- dimensional 
perturbed dynamical systems of a Newtonian type 

± = y, 

dV _ (12) 

y = --g^ix) +ap"'{x,y)y. 

We have used Belinskii parameterization of the viscosity 3£,{p) — ap™ and the first integral p — = ik/x"^ or 
Pcff — is preserved also in the case of viscous cosmology. If a = we obtain the standard conservative FRW 
cosmology. 

It is useful to distinguish two special cases. 



A. Flat models (k=0) 

In this case the dynamical effects of viscosity are equivalent to the effects of Chaplygin gas [H, [s^. If we put 
PcS = IP ~ Sap^^'H then from the conservation condition p = —3H{p + Pcs) we obtain the relation 

= (tT^ + :,3(i+f)(i+a) ) ' (13) 

where A and B are integration constants, and p{x = 1) = po = + B. 

In the special case of a = (m — —1/2) and 7 = (dust) the corresponding equation relating energy density to 
the scale factor for Chaplygin gas can be recovered. Therefore all the flat models with viscosity following Belinskii 
parameterization can by represented as conservative systems. To construct the phase space portrait it is useful to 
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represent the potential function V{x) in the similar way as for generalized Chaplygin gas for which the equation of 
state assumes the following form 



P 

and energy density can be parameterized by po and As 



A_ 



^(-) - PO + ^3(t,)f,;,) j , (14) 

where — aAg is the squared velocity of sound. Therefore the potential function for the system is given by 

1 

= -^- -if^viscOX^ + ^3(l;)fi;.) ) ' (15) 

where 1 + a = 1/2 - m, A^ = , f^visc.o = gfr- 

In Figs. [7| and [5] the phase portraits for different parameter m and for different equation of state coefficient 7 are 
illustrated. Therefore the effect of viscosity in the perfect fluid is equivalent to the effect of Chaplygin gas in the case 
of flat model. 

Note that for phantom cosmology (7 < — 1) Ag can be either negative or positive but if Ag is negative, the expression 
in the parentheses should be positive. Therefore we obtain following inequalities which restrict 



irnax 



_ / < (t) if (1 + 7)(1 + a) > 0, 

I > (^)'"+"'"^°' if (l+7)(l + a) < 0. 

Therefore for < (1 + 7 < 0) x < x^ax if only a < — 1 (m > 1/2) and x > x^i^i if only a > — 1 (m < 1/2). 

If we assume that 7 > — 1, then the Chaplygin gas interpolates between the state of positive cosmological constant 
and state of matter domination p cx a~^(^+'''^ at the early epoch. Hence the viscosity effects can shift structurally 
unstable model into open dense set. Fig.[9]shows the influence of viscosity on bouncing models for which the potential 
function contains additional contribution which can unify matter (p = 0) and dark energy. 

If we postulate the existence of negative contribution to the effective energy density coming, for example, from 
quantum effects like in Vandersloot's approach [13, HI], then the potential function assumes the following form 




Vix) = --<^ rivisc.ox^ As + —j-^ - A^ox-^ . (16) 



where r2visc,o ~ ^iq.o = 1 and then the cosmological model exhibits characteristic bounce phase [39|, |40|. 

It is characteristic for loop quantum cosmology that initial singularity is replaced with bounce because of quantum 
effects which are manifested by presence of negative energy contribution to the H'^{a) relation. In this picture space- 
time does not end at singularities and quantum regime can offer a bridge between vast space-time regions which are 
classically unrelated (see ^27)]). 

Unfortunately, in the phase space we obtain a non-hyperbolic critical point around the bounce which makes the 
system structurally unstable. In Fig.[n]it is illustrated how small viscosity effects (or equivalently the Chaplygin gas) 
can dramatically change global phase portraits and makes the corresponding system (with potential pop ) structurally 
stable. 

B. Case of constant viscosity (ni=0) 

Let us start now the analysis of the effects of viscosity in the case of constant viscosity. In this special case the 
equation of motion assumes a very simple form 

X ~ y, 

■ (^^) 
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c) 

FIG. 7: The phase portraits for the viscous cosmological models, a) 7 = 0, a = (m = —1/2), As = 0.3, b) 7 = 1/3, a = 1/2 
(m — —1), As = 0.5, c) 7 = —4/3, a = {m = —1/2), As = —0.95. Note that both models presented in Fig. a) and b) are 
topologically equivalent. 

One can calculate the distance from a dissipative system P7|) to a conservative one. Then the perturbation vector is 
and distance 

d{cons,dissip) = max {|ay|,Q;|=Q! max {|?;|,l}, 

{x,y)£E {x,y)£E 

where E is a compact region of the phase space in which we compare models, which we chose as the region \y\ < 1. 
Hence the parameter a will measure the distance from a viscous perturbed system to a conservative system, in 




FIG. 8: The phase portraits for the viscous cosmological models, d) 7 = —4/3, a — 1/2 (m = —1), As = —0.5, and the solid 
(blue in el. version) line denotes the limit x = Xmin, e) 7 = —4/3, a = —3/2 (m = 1), As — —0.5, the solid line denotes the 
hmit X = a;max. 



metric. It is interesting that the exact solution of the system (jl7p with constant viscosity can be simply obtained, 
vi^ithout any assumption about the curvature, if we know the corresponding solution without the viscosity. The 
following theorem establishes this fact. 

Theorem 1 Let t = 4i{a) he a solution of the conservative system with dark energy 
It can be determined from the relation 

dx 



tia) 



Then solution of the type 



exp [at] — (f){a), 



0(a). 



(18) 



will be the solution with constant viscosity a. 

To prove this fact, it is sufficient to check that solution of the above type will satisfy equation 

dV 



da 



aa, 



which is equivalent to equation p7p . 

One can consider now the properties of the phase plane of the dissipative dark energy models. We are especially 
interested whether bouncing models can be shifted to the generic set of the multiverse. 

The location of all critical points at finite domain coincides with that of the case without dissipation effects because 
all should be situated on the a;-axis (static critical point). Their character is determined by the linearization matrix 



A = 



1 



(19) 



(0:0,0) 
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FIG. 9: The bouncing models with viscosity for potential (|16|) with As — 0.75. f2visc,o > 1 (left) and fivisc.o < 1 (right). Note 
that viscosity makes the global phase portraits of bouncing models structurally stable. 

Because tr A = a > all the critical points are unstable. The characteristic equation assumes the following form 

\^-aX + V^,ixa) = 0. (20) 

Therefore the eigenvalues of the linearization matrix (I19p are real of the same signs if V^xixo) > and opposite signs 
if A1A2 = Vxx{xo) < 0. They are real if A = — AVxxixo) is positive, or imaginary if A is negative. In any case 
ReA = a is different from zero, i.e., all critical points are hyperbolic. Note that if Vxx{xo) < then always A > 0, 
i.e., eigenvalues are real of opposite signs (they are representing saddles). If Vxx{xo) > (like in the neighborhood of 
bounce) then eigenvalues are real and positive (Ai + X2 = a and A1A2 > 0) or imaginary (if Vxxixo) > a^/4). This 
corresponds to the presence of an unstable node or an unstable focus on phase portraits, respectively. Because both 
unstable nodes are topologically equivalent to the unstable focus, we obtain generic phase portraits for the dissipative 
bouncing model which are demonstrated in Fig. 1101 

Let us now present some general properties of the dissipative FRW models. The equation of motion is in the form 

i = 2/ = Pix,y), 

dV (21) 

y = ^"^(^) + "2/ = Q{x,y)- 

The critical points at finite domain (only repellors) are determined from 
System ((2T|) linearized around the static critical point (a:o,0) has the form 

{x - xo) + ( X" ) ^ 2^")' 

(:^^o,0) J {x„fi) 

{.x-xo)^[ — \ (y-yo), 

{xo,0} \Oy J (xofi) 



' ix 






- a;o)' = 


(-) 

\dx J 


< 

{y 


- 2/0 )■ = 


fdQ\ 
[dxJ 



dx'^ 



{X - Xq). 
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FIG. 10: The phase portrait of the system psp for the bouncing ACDM model with a constant viscosity coefficient which 
assumes different values. Cases: a) flA,o ~ 0.4, r2m,o = 0.3, a = 0.1, b) the same but with a = 0.9; c) JIa.o ~ 0.8, f2m,o ~ 0.3, 
a = 0.1, d) the same as previous but with a = 0.9. Note that all phase portraits are structurally stable due to viscosity 
perturbation and that symmetry x — > —x can be broken if viscosity effects are included. It is because entropy is changing due 
to the presence of dissipative processes. 



It is useful to shift the critical point to the origin of a new coordinate system labelled as {X, Y), say x ^ X = x — xq, 
Y = y. Then 
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or X = AX. Let us consider the special case of a conservative system for a = 0. The solutions of the system are 

X ^Ci exp [Xit] + C2 exp [\2t], 

where Ci and C2 are constants and Ai, A2 are eigenvalues of the linearization matrix. 
There are two types of the solutions: 



1. solution of a saddle point if eigenvalues are real of opposite signs A1.2 = ±\/— det A 



2. solution of a center type when eigenvalues are imaginary Ai_2 = det A — ±i 

Finally we obtain corresponding solution 



. - .0 = exp J ( - t] + C2 exp [J ( - t], 



x-xo = CiCos^(g) t + C,sinJ(^^ t. 

y \ / a;— a;o V x—xq 

In the case of dissipative [a = const) dynamical system the characteristic equation assumes the form ([^0]) and 
solutions are A1.2 — — "^^^ and A1A2 = det^ = {^^) x ' ^^^^ ^^'^ type of the critical points (if they exist) 
will depend on both the sign of discriminant A = — ^Vxx{xq) and sign of Vxx{xo) at the critical points, namely 



if Vxxixo) < 0, i.e. V{x) has a maximum then critical points are always saddles, 



• if Vxxixn) > 0, i.e. V{x) has a minimum at xq then character of the critical points will depend on the value of 
a such that: 



— if a > 2y^Vxx{xo) then we have unstable nodes, 

— if < a < 2^Vxx{xo) then we obtain unstable focus. 

Therefore in a is sufficiently small then unstable focus is present on the phase portrait instead of center in the 
conservative case. If a is sufficiently large then we obtain unstable node at the minimum of the diagram of the 
potential function. Note that except of degeneration case of center all admissible critical points are structurally 
stable. The focus is topologically equivalent to node. 

Note that the presence of dissipative effects enlarge the accelerating domain occupied by accelerating cosmological 
models because of relation 

dV 

-^ix) + ay>0, (22) 



and now this region depends additionally on the coordinate y. Note that at the critical point relation (|22p is always 
valid (like in [37|). 

If we substitute in (|17p the form of the potential function for bouncing cosmology with 7 = (dust) and n — A, 
then the corresponding dynamical system has the form 



x^y, 

1 -2 -3 (23) 

y^--^m,ox +^n,OX +flA,QX + ay. 

Rescaling the time variable we can make r.h.s. of the system (j23p polynomial. The simple time transformation 
given by T I— > 77: dr/x^ = drj give rise to its regularization. Such a choice of new time variable does not change 
the phase portrait because our system has been autonomous since the very beginning. This operation is equivalent 
to multiplication of the r.h.s. by x^ so that the system can be considered in the i] time, because x > and 7y is a 
monotonous function of any time t. Finally we obtain 

dX n 

- QjnflX + ^n,o + ^A,ox + ayx^ , 



I d?7 2 
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where n^.Q — ^n,o + ^A,o = 1- 

Note that reparameterization of cosmological time is equivalent to the choice of the corresponding lapse function 
N{t) in such a way that t — J N{t)dt. If we consider particle like representation of the dynamics then we obtain 
classical mechanics with lapse in the H.-J. Schmidt terminology [4]| . 

The degeneration of the dynamics at infinity for phantoms can be removed in a simple way by introducing a new 
positional variable x such that 



3 

(in general x — a;^(i+3™)/2) ^j^fj reparameterizing the time variable 

T f-* r/ : ^x^dt = ar], 
or : ^x^dt — drj. 

Of course the new time variable is a monotonous function of the original time r. Hence we obtain the new system in 
which the original variable y preserves its original sense 

da; 

d?7 ^' 

dy 1 / ^ 5 _N dV 

^ = --^ (^-iim.oa; ^ + SUphX J = ~'q^- 

New potential function assumes the form 

The above system is of Newtonian type and has a first integral in the energy conservation form 

^ + V{x) = 0, 

where now y = 

It is easy to check that the kinetic energy form = |(^)'^ is preserved under nonlinear rescaling transforma- 

tion X t—^ X and time reparameterization. The form of the potential function V{x) is de facto the same as the original 
V{x) in which instead of x we operate on x as the positional variable. 

From the first integral we obtain that for large x we have 



because material contribution is negligible. It means that asymptotically the trajectories reach the point which lies 
on the circle at infinity and on the straight line y oc x. The corresponding dynamical system on the plane (x, y) is 
of course equivalent to the ACDM model phase portrait although we must remember that the attractor at infinity 
does not represent the de Sitter universe. Therefore the obtained global phase portrait for phantom cosmology is 
structurally stable in a similar way as in the case of the ACDM model. 

The analogous "regularization procedure" which we perform for phantom with w = — | can be realized in a more 
general case for any w. It is sufficient to choose 

_ 1 + 3™ 

x I— > x : X — X 2 ^ 

T^ri: d?7 = -(ii|2i)2;-i(i+"''di. 



The corresponding potential function is 

V{x) = --|rim,oS"^ + ^w.Qx"^ + ^k,0 

Therefore phantom cosmology in which a weak energy condition is violated forms an open and dense subsets in the 
multiverse of all dark energy models. Also big-rip singularities which are attributed to these classes of models are a 
generic features of their long term behavior. 
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The universe is accelerating in the region of phase space determined by the condition 

- flmfiX + 2\ iln,o + ^A,ox^ + oyx^ \ > 0, (25) 



and finally region of accelerating expansion is given by 

Therefore all terms apart from the matter term fim.o act toward the acceleration of the Universe. 

IV. VISCOUS COSMOLOGY TESTED BY DISTANT SUPERNOVAE TYPE lA 

Using the energy conservation condition we can obtain a simple relation 

p(a) = 



A B 
+ 



1+7 a3(i+T)(5-™) 

where ^ = 3a is a positive constant which quantifies if viscous effects appeared, B is an arbitrary integration constant. 
For small value of the scale factor a{t) we have 

5l/(i-m) 

a3(i+7) ' 

which corresponds to the universe dominated by matter satisfying equation of state p = 7p. Also for a large value of 

scale factor (if 1 + 7 > 0) 



A 



and p ~ — I I = —p 



1 + 7/ Vl + 7. 

which corresponds to an otherwise empty universe with cosmological constant (j^)^^^^ For an accelerating 
universe, the deceleration parameter must be negative, therefore 

^sii+,)a-m) > i?(l + 7)(l + 37) 
3^ 

The above expression indicates that if (^ — m)(l + 7) > then a universe starts to accelerate when the scale factor 
is reaching its critical value 



B(l+7)(l + 37) 
2A 



3(l/2-m)(l+7) 



Hence the interpolation of early matter domination phase and cosmological constant epoch is achieved for viscous 
fluid when the viscosity decreases when density is decreasing. If we put a = 1 at the present epoch then present 
Universe is accelerating if (^)(1+7)(1+ 37) < 1. 

An important test allows to verify if the viscosity fluid may represent dark energy is the comparison with the 
supernovae type la data. For this purpose we have to calculate the luminosity distance in the model 

The Friedman equation can be rearranged to the form giving explicitly the Hubble function H{z) = a/a 



+ zf + rivisc + (1 - ^)(i + " 



(27) 



where the quantities fli, i = m, vise, k represent fractions of critical density currently contained in energy densities of 
respective components and f^m + f^visc + = 1- 
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This model adds to our previous discussion of a viscous fluid the additional elements of curvature and matter scaling 
like dust. As discussed above (see equation P^ ). the viscous fluid has a Chaplygin-like dependence on scale factor 
only in the flat dust-free case; in other circumstances, equation (j27p does not represent the solution for a universe 
containing in part a viscous fluid. Therefore, this model is in general to be regarded as an empirical construct. Such 
a model can nevertheless be compared to data, and we do this below; where we mention model testing, we actually 
mean the estimation of the m and As parameters for the best fit FRW cosmological model. We will mostly use the 
flat model A: = (the exception takes place when we relax a flat prior) since the evidence of this case is very strong 
in the light of current CMBR data. 

Formula (|28|) is the most general one in the framework of the Friedmann- Robertson- Walker models with this unified 
macroscopic (phenomenological) description of both dark energy and dark matter (quartessence) . We assume that 

Finally the luminosity distance reads 

di(z)-(l + z)— / , (28) 

^Jn,^.o{^ + zf + rivisco + (i - ^)(i + — 

To proceed with fitting the SNIa data we need the magnitude-redshift relation 

m(z, A^, fii„,o, f^visco; As,a) ^ M + b log^o Dl{z, fJm,o, f^visc,o; As,a) (29) 

where: 

Dl{z, fim,o, f^visc^o; As, a) = HodLiz, Hq, rim,o, ^Ch,o\ As, a) 

is the luminosity distance with Hq factored out so that marginalization over the intercept 

M-5 logio Ha + 25 (30) 

leads actually to joint marginalization over Hq and M [M being the absolute magnitude of SNIa). 
Then we can obtain the best fitted model minimizing the function 



v-^ (mf ' - mf^ 



,2 



where the sum is over the SNIa sample and Ui denote the (full) statistical error of magnitude determination. This is 
illustrated in Fig. [11] of residuals (with respect to the Einstein-de Sitter model) and levels in the {As,m) plane. 
One of the advantages of residual plots is that the intercept of the m — z curve gets cancelled. The assumption that 
the intercept is the same for different cosmological models is legitimate since M is actually determined from the 
low-redshift part of the Hubble diagram which should be linear in all realistic cosmologies. 

The best fit values alone are not relevant if not supplemented with the confidence levels for the parameters. 
Therefore, we performed the estimation of model parameters using the minimization procedure, based on the likelihood 
function. We assumed that supernovae measurements came with uncorrelated Gaussian errors and in this case the 
likelihood function C could be determined from chi-square statistic C oc exp (— x^/2) [J, 

Therefore we supplement our analysis with confidence intervals in the {As , m) plane by calculating the marginal 
probability density functions 

V{As,a)oc J exp(^~^x^{n^fi,n^isc,o,As,a,M)^dM 

with rJni^Oj f^visc,o fixed (flm.Q = 0.0,0.05,0.3) and 

V{As,a) cx y exp ^ - ^x^i^m.o, revise, o, ^s, a, M)^ dil^.o 

with Ai fixed {A4 — —3.39) respectively (a proportionality sign means equal up to the normalization constant) 
fFigll2p . In order to complete the picture we have also derived one-dimensional probability distribution functions for 
revise obtained from joint marginalization over m and As- The maximum value of such PDF informs us about the 
most probable value of f2visc,o (supported by supernovae data) within the full class of viscous cosmological models. 
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FIG. 11; Results from the Gold sample. Residuals (in mag) between the Einstein-de Sitter model (zero line), and three flat 
models: with m — 0, the flat ACDM model (m — —0.5), and m = —1.5 (equivalent to Chaplygin gas model with a = 1). Note 
that in the viscous cosmology distant z > 1 supernovae should be brighter than in the ACDM model. Therefore future distant 
supernovae data can definitely distinguish; which cosmology, with dark energy with dissipation or not? 



Adopting an analogy to the models with generalized Chaplygin gas (GCG) one can see that pvisc,o — (1^ + -^) " 
represents the current energy density of the viscous fluid for which one can define density parameters in the standard 
way. One can also calculate the adiabatic speed of sound squared for viscous fluid 

2 _ 2 '' 1+7 

'^s — A I B ■ 

1+7 ~^ a3<i/2-m) 

Therefore q = " '"^^ (velocity of light at present epoch) and constants A and B can be expressed in terms of 

quantities having well defined physical meaning. Instead of constants A and B it is useful to consider new ones 
As = j^, B = 1 — As and pvisc.o- Then cj. Q = (— 5 — m)As in the units of the speed of light c. 



A. Fits to As and m parameters 



1. Samples used 

We use data sets which were compiled by Riess et al. [l^ , they have been used by many researchers as a standard 
dataset. They improved the former Riess et al. sample and discovered 16 new type la supernovae. It should be noted 
that 6 of these objects have z > 1.25 (out of total number of 7 object with so high red shifts). Moreover, they compiled 
a set of previously observed SNIa relying on large, published samples, whenever possible, to reduce systematic errors 
from differences in calibrations. Thanks to this enriched sample it became possible to test our prediction that distant 
supernovae should be brighter in viscous cosmology than in the ACDM model (see discussion below) . 

The full Riess sample contains 186 SNIa ("Silver" sample). Taking into account the quality of the spectroscopic and 
photometric record for individual supernovae, they also selected a more restricted "Gold" sample of 157 supernovae. 
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FIG. 12; Confidence levels on the (j4s,m) plane for the viscous cosmological model, marginalized over A4, and fim,o- 

We have separately analyzed the ACDM model for supernovae with z < 1 and for all SNIa belonging to the Gold 
sample. 

2. Viscous cosmological models tested 

Using these samples we have tested viscous cosmology in three different classes of models with (1) ilmfl — 0.3, 
^^visc,o = 0.7; (2) r2m,o = 0.05, r2visc,o — 0.95 and (3) r^m.o — 0, revise, o — 1- We started with a fixed value of 
M — —3.39 modifying this assumption accordingly while analyzing different samples. 

The first class was chosen as representative of the standard knowledge of r^m.o (baryonic plus dark matter in galactic 
halos with viscosity responsible for the missing part of closure density (the dark energy). 

In the second class we have incorporated (at the level of rim.o) the prior knowledge about the baryonic content of 
the Universe (as inferred from the BBN considerations). Hence this class is representative of the models in which 
the viscous fluid is allowed to clump and is responsible both for dark matter in halos as well as its diffuse part (dark 
energy). 
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The third class is a kind of toy model - the FRW universe filled completely with viscous fluid. We have considered 
it mainly in order to see how sensitive the SNIa test is with respect to parameters identifying the cosmological model. 
Finally, we analyzed the data without any prior assumption about fimjO- 

3. Results 

For statistical analysis we restricted the values of the As parameter [0,1.0] and m to the interval [—0.5,-1.5] 
because of the relation = aAg (where a — —0.5 — to). This constraint guarantee that in all cases c^. is real and 
does not exceed c = 1. We separately analyze the case of m = 0, which corresponds to constant viscosity coefficient. 

The results (best fits) of two fitting procedures performed on Riess samples and with different prior assumptions 
concerning the cosmological models are presented in Tables |TT] and IIIII Table |lT] refers to the method whereas 
in Table IIIII the results from marginalized probability density functions is displayed. In both cases we obtained 
different values of M for each analyzed sample. In Table HV] we gather results of the statistical analysis of the viscous 
cosmological model performed on the Gold Riess sample of SNIa as a minimum best fit for different values of fixed 

TO. 

First, we present residual plots of redshift-magnitude relations between the Einstcin-de Sitter model (represented 
by zero line) and models with to = — upper curve, m — —0.5 (this model is equivalent to the ACDM model) — 
middle curve, and m — —1.5 (this model is equivalent to the GCG) — lower curve. One can observe that systematic 
deviation between these models is larger at higher red-shifts. The viscous model (m = —1.5) predicts that high 
redshift supernovae should be brighter than what is predicted with the ACDM model, while in the model with m = 
high redshift supernovae should be fainter than those predicted with the ACDM model. 

The Riess sample leads to the results which are similar to these obtained recently with the Astier sample [44| . 
For the Gold sample, the joint marginalization over parameters gives the following results: ilvisc,o = 1-00 (hence 
^m,Q = 0.0), with the limit rJviscO > 0.80 at the confidence level of 68.3% and fivisco > 0.69 at the confidence level 
of 95.4%. (to = -1.5, As = 0.83) with the interval to e (-0.86,-1.5) and As £ (0.76, 0.94) at the confidence level of 
68.3% and to e (-0.55, -1.5) and As S (0.72, 1.00) at the confidence level of 95.4%. 

B. Information criteria 

In the modern observational cosmology the so called degeneracy problem is present: many models with dramatically 
different scenarios are in good agreement with the present observations. Information criteria of the model selection 
[45| can be used to solve this degeneracy among some subclass of dark energy models. Among these criteria the 
Akaike information criterion (AIC) jl^ and the Bayesian information criterion (BIG) [47j are most popular. From 
these criteria we can determine the number of the essential model parameters providing the preferred fit to the data. 

The AIC is defined in the following way [46l ] 

AIC = -21n/: + 2d (31) 

where C is the maximum likelihood and d is a number of the model parameters. The best model with a parameter 
set providing the preferred fit to the data is the one that minimizes the AIC. 
The BIG introduced by Schwarz [l^ is defined as 

BIG = -21n£ + dlniV (32) 

where A'' is the number of data points used in the fit. While the AIC tends to favor models with large number of 
parameters, the BIG penalizes them more strongh^, so the BIG provides a useful approximation to full evidence in the 
case of no prior on the set of model parameters [4^] . 

The effectiveness of using these criteria in the current cosmological applications has been recently demonstrated 
by Liddle [1^ who, taking GMB WMAP data ^], found the number of essential cosmological parameters to be five. 
Moreover he obtained the important conclusion that spatially-flat models are statistically preferred to close models 
as it was indicated by the GMB WMAP analysis (their best-fit value is il.tot,o = '^i^i,Q — 1-02 ± 0.02 at la level). 

In the paper of Parkinson et al. ^SJ the usefulness of Bayesian model selection criteria in the context of testing 
for double inflation with WMAP was demonstrated. These criteria was also used recently by us to show that models 
with the big-bang scenario are rather preferred over the bouncing scenario [39j . 

Please note that both information criteria values have no absolute sense and only the relative values for different 
models are statistically interesting. For the BIG a difference of 2 is treated as a positive evidence (6 as a strong 
evidence) against the model with larger value of the BIG [13, HH . Therefore one can order all models which belong 
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TABLE II: Results of the statistical analysis of the viscous cosmological model (with marginalization over A^) performed on 
the Riess samples of SNIa (Silver, Gold) as a minimum ^ best fit. First rows for each sample refer to no prior on Jlm.o- The 
same analysis was repeated with fixed priors Slm,o ~ 0.0, fim.o = 0.05 and r^m.o — 0.3. 



sample 




i ^visc,0 


A 

As 


m 




2 

X 


Silver 


0.00 


1.00 


0.82 


-1.50 


15.945 


229.4 




0.00 


1.00 


0.82 


-1.50 


15.945 


229.4 




0.05 


0.95 


0.85 


-1.50 


15.945 


229.6 




0.30 


0.70 


0.99 


-1.50 


15.965 


232.3 


Gold 


0.00 


1.00 


0.81 


-1.50 


15.945 


173.7 




0.00 


1.00 


0.81 


-1.50 


15.945 


173.7 




0.05 


0.95 


0.84 


-1.50 


15.945 


173.8 




0.30 


0.70 


0.99 


-1.50 


15.965 


175.6 



TABLE III: viscous cosmological model parameter values obtained from the marginal probability density functions calculated 
on the Riess samples. First rows for each sample refer to no prior on f2m,o- The same analysis was repeated with fixed priors 
f^m.o = 0.0, Qmfi = 0.05 and fim.o ~ 0.3. 



sample 






As 


m 


M 


Silver 


0.00+'''^« 


1.00-0.18 


U-04_o 06 


_150.59 


i5.945j:E;:| 




0.00 


1.00 


n 70+0.03 

'y-0.05 


-1.5°-S2 


15.955 Q 




0.05 


0.95 


n 01+0.04 


-1.50-5* 


15.955 Q Q2 




0.30 


0.70 


n qq+0.01 

U.»»_o 03 


— 0.5-0.64 


15.965tH2 


Gold 


O.OO+^-^o 


1.00-0.20 


0.83lo;j7 


_l 5+0.64 


15.955 Q Q2 




0.00 


1.00 


0.77lH^ 


-l-5o.58 


15.955 Q Q2 




0.05 


0.95 


o.soint 


-l-5o.59 


15.955 Q 




0.30 


0.70 


0.99lHi 


_Q5-0.64 


15.9651H2 



to the ensemble of dark energy models following the AIC and BIC values. If we do not find any positive evidence 
from information criteria the models are treated as identical and eventually additional parameters are treated as not 
statistically significant. Therefore the information criteria offer the possibility of introducing relation of weak order 
in the considered class of analyzed models. 

In Table IVl we present results of analysis of four flat dark energy models with two free parameters, i.e., the models 
with m — 0, —0.5, —1 — 1.5. On could note that both the AIC and BIC prefer the case m — —1.5 i.e viscosity model 
with m — —1.5 over the ACDM model. It is interesting to observe that if the value of m parameter was derived from 
physics then number of independent model parameter would be lower by one and the viscous model is favored over 
the ACDM model by model selection criteria (the values of AIC and BIC for the flat ACDM model are 179.9 and 
186.0, respectively, see j52i]). 

1. Bayes factor 

For deeper statistical analysis it is used the Bayes factor [H, [sj] ■ Using this technique we compare the viscous 
models with the concordance ACDM model. 



TABLE IV: Results of the statistical analysis of the viscous cosmological model (with marginalization over A4) performed on 
the Gold Riess sample of SNIa as a minimum best fit for different values of fixed m. 



sample 


f2m,0 


^visc,0 


As 


m 


M 


x' 


Gold 


0.05 


0.95 


0.65 


0.0 


15.965 


177.9 




0.05 


0.95 


0.73 


-0.5 


15.955 


175.9 




0.05 


0.95 


0.79 


-1.0 


15.955 


174.6 




0.05 


0.95 


0.84 


-1.5 


15.945 


173.8 
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TABLE V: Results of the AIC and BIC for models with different values of m (with marginalization over A4) performed on Gold 
Riess sample of SNIa. Note that if the m parameter is fixed by physics of the viscous process then viscous model is preferred 
over the ACDM one. 

sample m AIC BIC 
Gold 0.0 183.9 188.0 

-0.5 179.9 186.0 

-1.0 178.6 184.7 

-1.5 177.9 183.9 



In the Bayesian framework to compare models (the model set {Mi}, i — 1, . . . , K) is to find the value of probability 
in the light of data (so called a posterior probability) for each model. We can define the posterior odds for models 
Mi and Mj , which (in the case when no model is favored a priori) is reduced to the marginal likelihood {E) ratio (so 
called the Bayes factor - Bij) 

_ jL{e\D,M,)pmu)de ^ 

jmD,M,)P{fj\M,)dv £;/ ^ > 

where 6 is the parameter vector, which defines model i, L{6\D, Mi) is the likelihood under a model i, P{9\Mi) is the 
prior probability for 6 under a model i. 

It is interesting to compare the flat FRW model with bulk viscosity and fixed baryonic matter {ilm.o — ^h,a — 0.05) 
with the concordance ACDM model. 

To compare models Mi and Mj one can compute 2\nBij — — (BIC^ — BICj) = — ABICy which can be interpret 
as 'strength of evidence' against j model: < 2 In By < 2-not worth more than a bare mention, 2 < 2 In Bij < 6 - 
positive, 6 < 2\n Bij < 10 - strong, and 2ln Bij > 10 - very strong. 

It is useful to choose one model from our model set (a reference model-s) and compare the rest models with this 
one model, situation in which 2hiBsi > indicates evidence against model i with respect to the reference model, 
whereas 2 In B si < denotes evidence in favor of model i. 

We can compute posterior probability for model i in the following way 

P{M.\D) = , (34) 

l^k=l Bks 

where Bis = exp[iAi3/Csi]. Then one can see how prior believe about model probability P{Mi) = change after 
inclusion data to analysis. This is the probability for model i being the best model from the set of models under 
consideration. 

Assuming equal priors for both concurrence models 0.5 we can calculate the posterior probabilities as 0.7408 for 
the viscous model and 0.2592 for the ACDM model. This indicates that the former model should be treated seriously 
as a candidate for dark energy description (see also [36l|). 



V. CONCLUSIONS 



In this paper we apply the theory of qualitative investigations of differential equations to the study of dissipative 
cosmological model with dark energy. We show that application of qualitative theory of dynamical system allows to 
reveal some structural stability properties of this model. 

We have developed a phenomenological unified model for dark energy and dark matter through the dissipation 
effects acting in the flat FRW model. 

We demonstrate that they describe a smooth transition from a decelerated expansion phase dominated by matter 
contribution to the present dark energy epoch. We can find one-to-one correspondence between models with general- 
ized Chaplygin gas {p = —A/p") and models filled by viscous fluid (p = jp — 3S,{p)H) in the flat case. It is established 
if we use Belinskii power law parameterization of viscosity coefficient ^(p) cx and 1 + a — ^ — m. 

We conclude that while the ACDM model comes in good agreement with other observational data, mainly those 
from WMAP and dynamics of clusters of galaxies, the fitting quality for dissipative model is comparable to that value 
obtained for models with A. 

We supplemented our analysis with confidence intervals in the {As,m) plane. Our result show also that viscous 
cosmology predicts that at distance z > 1 supernovae should be brighter than in the ACDM model. 
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For deeper analysis of statistical results and to decide which model (with dark energy with or without dissipation) is 
distinguished we use the Akaike and Bayesian information criteria. Applying the model selection criteria we show that 
both the AIC and BIC indicate that additional contribution arising from the dissipative effects should be incorporated 
to the model if only physics fix the value of the parameter m. 

Since classical papers of Murphy [1^] and Heller et al. [13] it is well known that bulk viscosity can produce 
cosmological models without the initial singularity for flat universes. In an analogous way one can consider influence 
of bulk viscosity effects on avoiding the singularities in the future, namely so called big-rip singularities appearing in 
phantom cosmology [H^. Let us note that if A^. < (7 < — 1 and A > 0) then we have upper bound on the value 
of the scale factor if only m > 1/2. This means that viscosity effects acting in phantom's matter can give rise to 
avoidance of big-rip singularity in the future. 

Using the analogy of dissipative cosmology to the conservative one with the Chaplygin gas, one can find at least 
three significant features of dissipative cosmology |56j. First, similarly to the conservative models with the Chaplygin 
gas they describe a smooth transition from a decelerated to the present accelerated expansion. Second, the models 
attempt to give a unified phenomenological description of both dark energy and dark matter in a natural way without 
references to effects of immersion of our universe into higher dimensional bulk spaces or tachyon cosmological models. 
Third, they are a simple and natural extension of the ACDM model which has the property of flexibility with regard 
to the observational data. 

The main conclusions from the dynamical systems analysis are the following 

• The general FRW model with bulk viscosity parameterized by the Belinskii power-law parameterization is the 
two-dimensional dynamical system and framework of dynamical systems methods can be applied; we showed 
that the phase portrait of the model with bulk viscosity with parameter values obtained from statistical analysis 
is equivalent to the phase portrait of the ACDM model. 

• For the special case of the flat model the effect of bulk viscosity are formally equivalent to the effects of 
generalized Chaplygin gas. We have found simple relation between parameter a from the equation of state of 
the generalized Chaplygin gas {p — —Ap^") and the parameter m taken from the Belinskii parameterization 
of viscosity ^ oc p™. These models are formally equivalent and cannot be distinguished using the kinematic 
cosmological tests based on the observables obtained from null geodesies. 

• Moreover in the case of flat models dynamics can be reduced to the form of conservative system of the Newtonian 
type. 

• In the special case of the flat model with constant viscosity we have formulated the theorem 1 which states that 
the solution for viscous model can be obtain from the the solution of the model without the viscosity. 

An additional argument for incorporating of the new dissipative parameter is brought by the theory which will always 
favor the structurally stable models over fragile ones. We show that viscosity effects give rise to structurally stable 
evolutional scenario with squeezing bounce phase predicted by loop quantum gravity. However we must remember 
that flat cosmological model with viscosity effects is equivalent to the FRW model with the generalized Chaplygin 
gas and is not better than the ACDM model. 

There is a class of dynamical systems of cosmological origin like the cosmological models of loop quantum gravity 
in which instead of an initial singularity we have a bounce [57| . The corresponding dynamical systems is structurally 
unstable (see Section II) because of the presence of the center in the phase portrait. Therefore following the Peixoto 
theorem they are exceptional among the cosmological models on the plane. In the bouncing cosmology the acceleration 
of the Universe is only transient phenomenon and this is a non-generic case. Due to effect of bulk viscosity they become 
structurally stable models. 

One can derive some philosophical conclusions from our analysis. The arrow of time [58] is a physical mystery 
because while fundamental laws of physics are CPT invariant there is thermodynamical arrow of time prescribed by 
the second law of thermodynamics and cosmological arrow of time prescribed by expansion of the Universe. If our 
Universe is not flat then dissipative effects of bulk viscosity can determine the cosmological arrow of time. 
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